home *** CD-ROM | disk | FTP | other *** search
/ Australian Personal Computer 2002 November / CD 1 / APC0211D1.ISO / workshop / prog / files / ActivePerl-5.6.1.633-MSWin32.msi / _fb3d694d9d929bd5adcfcc4517b4f5f6 < prev    next >
Encoding:
Text File  |  2001-09-04  |  14.1 KB  |  520 lines
  1. package Math::BigInt;
  2. $VERSION='0.01';
  3.  
  4. use overload
  5. '+'    =>    sub {new Math::BigInt &badd},
  6. '-'    =>    sub {new Math::BigInt
  7.                $_[2]? bsub($_[1],${$_[0]}) : bsub(${$_[0]},$_[1])},
  8. '<=>'    =>    sub {$_[2]? bcmp($_[1],${$_[0]}) : bcmp(${$_[0]},$_[1])},
  9. 'cmp'    =>    sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
  10. '*'    =>    sub {new Math::BigInt &bmul},
  11. '/'    =>    sub {new Math::BigInt 
  12.                $_[2]? scalar bdiv($_[1],${$_[0]}) :
  13.              scalar bdiv(${$_[0]},$_[1])},
  14. '%'    =>    sub {new Math::BigInt
  15.                $_[2]? bmod($_[1],${$_[0]}) : bmod(${$_[0]},$_[1])},
  16. '**'    =>    sub {new Math::BigInt
  17.                $_[2]? bpow($_[1],${$_[0]}) : bpow(${$_[0]},$_[1])},
  18. 'neg'    =>    sub {new Math::BigInt &bneg},
  19. 'abs'    =>    sub {new Math::BigInt &babs},
  20. '<<'    =>    sub {new Math::BigInt
  21.                $_[2]? blsft($_[1],${$_[0]}) : blsft(${$_[0]},$_[1])},
  22. '>>'    =>    sub {new Math::BigInt
  23.                $_[2]? brsft($_[1],${$_[0]}) : brsft(${$_[0]},$_[1])},
  24. '&'    =>    sub {new Math::BigInt &band},
  25. '|'    =>    sub {new Math::BigInt &bior},
  26. '^'    =>    sub {new Math::BigInt &bxor},
  27. '~'    =>    sub {new Math::BigInt &bnot},
  28.  
  29. qw(
  30. ""    stringify
  31. 0+    numify)            # Order of arguments unsignificant
  32. ;
  33.  
  34. $NaNOK=1;
  35.  
  36. sub new {
  37.   my($class) = shift;
  38.   my($foo) = bnorm(shift);
  39.   die "Not a number initialized to Math::BigInt" if !$NaNOK && $foo eq "NaN";
  40.   bless \$foo, $class;
  41. }
  42. sub stringify { "${$_[0]}" }
  43. sub numify { 0 + "${$_[0]}" }    # Not needed, additional overhead
  44.                 # comparing to direct compilation based on
  45.                 # stringify
  46. sub import {
  47.   shift;
  48.   return unless @_;
  49.   die "unknown import: @_" unless @_ == 1 and $_[0] eq ':constant';
  50.   overload::constant integer => sub {Math::BigInt->new(shift)};
  51. }
  52.  
  53. $zero = 0;
  54.  
  55. # overcome a floating point problem on certain osnames (posix-bc, os390)
  56. BEGIN {
  57.     my $x = 100000.0;
  58.     my $use_mult = int($x*1e-5)*1e5 == $x ? 1 : 0;
  59. }
  60.  
  61. # normalize string form of number.   Strip leading zeros.  Strip any
  62. #   white space and add a sign, if missing.
  63. # Strings that are not numbers result the value 'NaN'.
  64.  
  65. sub bnorm { #(num_str) return num_str
  66.     local($_) = @_;
  67.     s/\s+//g;                           # strip white space
  68.     if (s/^([+-]?)0*(\d+)$/$1$2/) {     # test if number
  69.     substr($_,$[,0) = '+' unless $1; # Add missing sign
  70.     s/^-0/+0/;
  71.     $_;
  72.     } else {
  73.     'NaN';
  74.     }
  75. }
  76.  
  77. # Convert a number from string format to internal base 100000 format.
  78. #   Assumes normalized value as input.
  79. sub internal { #(num_str) return int_num_array
  80.     local($d) = @_;
  81.     ($is,$il) = (substr($d,$[,1),length($d)-2);
  82.     substr($d,$[,1) = '';
  83.     ($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d)));
  84. }
  85.  
  86. # Convert a number from internal base 100000 format to string format.
  87. #   This routine scribbles all over input array.
  88. sub external { #(int_num_array) return num_str
  89.     $es = shift;
  90.     grep($_ > 9999 || ($_ = substr('0000'.$_,-5)), @_);   # zero pad
  91.     &bnorm(join('', $es, reverse(@_)));    # reverse concat and normalize
  92. }
  93.  
  94. # Negate input value.
  95. sub bneg { #(num_str) return num_str
  96.     local($_) = &bnorm(@_);
  97.     return $_ if $_ eq '+0' or $_ eq 'NaN';
  98.     vec($_,0,8) ^= ord('+') ^ ord('-');
  99.     $_;
  100. }
  101.  
  102. # Returns the absolute value of the input.
  103. sub babs { #(num_str) return num_str
  104.     &abs(&bnorm(@_));
  105. }
  106.  
  107. sub abs { # post-normalized abs for internal use
  108.     local($_) = @_;
  109.     s/^-/+/;
  110.     $_;
  111. }
  112.  
  113. # Compares 2 values.  Returns one of undef, <0, =0, >0. (suitable for sort)
  114. sub bcmp { #(num_str, num_str) return cond_code
  115.     local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
  116.     if ($x eq 'NaN') {
  117.     undef;
  118.     } elsif ($y eq 'NaN') {
  119.     undef;
  120.     } else {
  121.     &cmp($x,$y) <=> 0;
  122.     }
  123. }
  124.  
  125. sub cmp { # post-normalized compare for internal use
  126.     local($cx, $cy) = @_;
  127.     
  128.     return 0 if ($cx eq $cy);
  129.  
  130.     local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1));
  131.     local($ld);
  132.  
  133.     if ($sx eq '+') {
  134.       return  1 if ($sy eq '-' || $cy eq '+0');
  135.       $ld = length($cx) - length($cy);
  136.       return $ld if ($ld);
  137.       return $cx cmp $cy;
  138.     } else { # $sx eq '-'
  139.       return -1 if ($sy eq '+');
  140.       $ld = length($cy) - length($cx);
  141.       return $ld if ($ld);
  142.       return $cy cmp $cx;
  143.     }
  144. }
  145.  
  146. sub badd { #(num_str, num_str) return num_str
  147.     local(*x, *y); ($x, $y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
  148.     if ($x eq 'NaN') {
  149.     'NaN';
  150.     } elsif ($y eq 'NaN') {
  151.     'NaN';
  152.     } else {
  153.     @x = &internal($x);             # convert to internal form
  154.     @y = &internal($y);
  155.     local($sx, $sy) = (shift @x, shift @y); # get signs
  156.     if ($sx eq $sy) {
  157.         &external($sx, &add(*x, *y)); # if same sign add
  158.     } else {
  159.         ($x, $y) = (&abs($x),&abs($y)); # make abs
  160.         if (&cmp($y,$x) > 0) {
  161.         &external($sy, &sub(*y, *x));
  162.         } else {
  163.         &external($sx, &sub(*x, *y));
  164.         }
  165.     }
  166.     }
  167. }
  168.  
  169. sub bsub { #(num_str, num_str) return num_str
  170.     &badd($_[$[],&bneg($_[$[+1]));    
  171. }
  172.  
  173. # GCD -- Euclids algorithm Knuth Vol 2 pg 296
  174. sub bgcd { #(num_str, num_str) return num_str
  175.     local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
  176.     if ($x eq 'NaN' || $y eq 'NaN') {
  177.     'NaN';
  178.     } else {
  179.     ($x, $y) = ($y,&bmod($x,$y)) while $y ne '+0';
  180.     $x;
  181.     }
  182. }
  183.  
  184. # routine to add two base 1e5 numbers
  185. #   stolen from Knuth Vol 2 Algorithm A pg 231
  186. #   there are separate routines to add and sub as per Kunth pg 233
  187. sub add { #(int_num_array, int_num_array) return int_num_array
  188.     local(*x, *y) = @_;
  189.     $car = 0;
  190.     for $x (@x) {
  191.     last unless @y || $car;
  192.     $x -= 1e5 if $car = (($x += (@y ? shift(@y) : 0) + $car) >= 1e5) ? 1 : 0;
  193.     }
  194.     for $y (@y) {
  195.     last unless $car;
  196.     $y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0;
  197.     }
  198.     (@x, @y, $car);
  199. }
  200.  
  201. # subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
  202. sub sub { #(int_num_array, int_num_array) return int_num_array
  203.     local(*sx, *sy) = @_;
  204.     $bar = 0;
  205.     for $sx (@sx) {
  206.     last unless @sy || $bar;
  207.     $sx += 1e5 if $bar = (($sx -= (@sy ? shift(@sy) : 0) + $bar) < 0);
  208.     }
  209.     @sx;
  210. }
  211.  
  212. # multiply two numbers -- stolen from Knuth Vol 2 pg 233
  213. sub bmul { #(num_str, num_str) return num_str
  214.     local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
  215.     if ($x eq 'NaN') {
  216.     'NaN';
  217.     } elsif ($y eq 'NaN') {
  218.     'NaN';
  219.     } else {
  220.     @x = &internal($x);
  221.     @y = &internal($y);
  222.     &external(&mul(*x,*y));
  223.     }
  224. }
  225.  
  226. # multiply two numbers in internal representation
  227. # destroys the arguments, supposes that two arguments are different
  228. sub mul { #(*int_num_array, *int_num_array) return int_num_array
  229.     local(*x, *y) = (shift, shift);
  230.     local($signr) = (shift @x ne shift @y) ? '-' : '+';
  231.     @prod = ();
  232.     for $x (@x) {
  233.       ($car, $cty) = (0, $[);
  234.       for $y (@y) {
  235.     $prod = $x * $y + ($prod[$cty] || 0) + $car;
  236.         if ($use_mult) {
  237.     $prod[$cty++] =
  238.       $prod - ($car = int($prod * 1e-5)) * 1e5;
  239.         }
  240.         else {
  241.     $prod[$cty++] =
  242.       $prod - ($car = int($prod / 1e5)) * 1e5;
  243.         }
  244.       }
  245.       $prod[$cty] += $car if $car;
  246.       $x = shift @prod;
  247.     }
  248.     ($signr, @x, @prod);
  249. }
  250.  
  251. # modulus
  252. sub bmod { #(num_str, num_str) return num_str
  253.     (&bdiv(@_))[$[+1];
  254. }
  255.  
  256. sub bdiv { #(dividend: num_str, divisor: num_str) return num_str
  257.     local (*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
  258.     return wantarray ? ('NaN','NaN') : 'NaN'
  259.     if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0');
  260.     return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0);
  261.     @x = &internal($x); @y = &internal($y);
  262.     $srem = $y[$[];
  263.     $sr = (shift @x ne shift @y) ? '-' : '+';
  264.     $car = $bar = $prd = 0;
  265.     if (($dd = int(1e5/($y[$#y]+1))) != 1) {
  266.     for $x (@x) {
  267.         $x = $x * $dd + $car;
  268.             if ($use_mult) {
  269.         $x -= ($car = int($x * 1e-5)) * 1e5;
  270.             }
  271.             else {
  272.         $x -= ($car = int($x / 1e5)) * 1e5;
  273.             }
  274.     }
  275.     push(@x, $car); $car = 0;
  276.     for $y (@y) {
  277.         $y = $y * $dd + $car;
  278.             if ($use_mult) {
  279.         $y -= ($car = int($y * 1e-5)) * 1e5;
  280.             }
  281.             else {
  282.         $y -= ($car = int($y / 1e5)) * 1e5;
  283.             }
  284.     }
  285.     }
  286.     else {
  287.     push(@x, 0);
  288.     }
  289.     @q = (); ($v2,$v1) = @y[-2,-1];
  290.     $v2 = 0 unless $v2;
  291.     while ($#x > $#y) {
  292.     ($u2,$u1,$u0) = @x[-3..-1];
  293.     $u2 = 0 unless $u2;
  294.     $q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1));
  295.     --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*1e5+$u2);
  296.     if ($q) {
  297.         ($car, $bar) = (0,0);
  298.         for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
  299.         $prd = $q * $y[$y] + $car;
  300.                 if ($use_mult) {
  301.         $prd -= ($car = int($prd * 1e-5)) * 1e5;
  302.                 }
  303.                 else {
  304.         $prd -= ($car = int($prd / 1e5)) * 1e5;
  305.                 }
  306.         $x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0));
  307.         }
  308.         if ($x[$#x] < $car + $bar) {
  309.         $car = 0; --$q;
  310.         for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
  311.             $x[$x] -= 1e5
  312.             if ($car = (($x[$x] += $y[$y] + $car) > 1e5));
  313.         }
  314.         }   
  315.     }
  316.     pop(@x); unshift(@q, $q);
  317.     }
  318.     if (wantarray) {
  319.     @d = ();
  320.     if ($dd != 1) {
  321.         $car = 0;
  322.         for $x (reverse @x) {
  323.         $prd = $car * 1e5 + $x;
  324.         $car = $prd - ($tmp = int($prd / $dd)) * $dd;
  325.         unshift(@d, $tmp);
  326.         }
  327.     }
  328.     else {
  329.         @d = @x;
  330.     }
  331.     (&external($sr, @q), &external($srem, @d, $zero));
  332.     } else {
  333.     &external($sr, @q);
  334.     }
  335. }
  336.  
  337. # compute power of two numbers -- stolen from Knuth Vol 2 pg 233
  338. sub bpow { #(num_str, num_str) return num_str
  339.     local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
  340.     if ($x eq 'NaN') {
  341.     'NaN';
  342.     } elsif ($y eq 'NaN') {
  343.     'NaN';
  344.     } elsif ($x eq '+1') {
  345.     '+1';
  346.     } elsif ($x eq '-1') {
  347.     &bmod($x,2) ? '-1': '+1';
  348.     } elsif ($y =~ /^-/) {
  349.     'NaN';
  350.     } elsif ($x eq '+0' && $y eq '+0') {
  351.     'NaN';
  352.     } else {
  353.     @x = &internal($x);
  354.     local(@pow2)=@x;
  355.     local(@pow)=&internal("+1");
  356.     local($y1,$res,@tmp1,@tmp2)=(1); # need tmp to send to mul
  357.     while ($y ne '+0') {
  358.       ($y,$res)=&bdiv($y,2);
  359.       if ($res ne '+0') {@tmp=@pow2; @pow=&mul(*pow,*tmp);}
  360.       if ($y ne '+0') {@tmp=@pow2;@pow2=&mul(*pow2,*tmp);}
  361.     }
  362.     &external(@pow);
  363.     }
  364. }
  365.  
  366. # compute x << y, y >= 0
  367. sub blsft { #(num_str, num_str) return num_str
  368.     &bmul($_[$[], &bpow(2, $_[$[+1]));
  369. }
  370.  
  371. # compute x >> y, y >= 0
  372. sub brsft { #(num_str, num_str) return num_str
  373.     &bdiv($_[$[], &bpow(2, $_[$[+1]));
  374. }
  375.  
  376. # compute x & y
  377. sub band { #(num_str, num_str) return num_str
  378.     local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
  379.     if ($x eq 'NaN' || $y eq 'NaN') {
  380.     'NaN';
  381.     } else {
  382.     while ($x ne '+0' && $y ne '+0') {
  383.         ($x, $xr) = &bdiv($x, 0x10000);
  384.         ($y, $yr) = &bdiv($y, 0x10000);
  385.         $r = &badd(&bmul(int $xr & $yr, $m), $r);
  386.         $m = &bmul($m, 0x10000);
  387.     }
  388.     $r;
  389.     }
  390. }
  391.  
  392. # compute x | y
  393. sub bior { #(num_str, num_str) return num_str
  394.     local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
  395.     if ($x eq 'NaN' || $y eq 'NaN') {
  396.     'NaN';
  397.     } else {
  398.     while ($x ne '+0' || $y ne '+0') {
  399.         ($x, $xr) = &bdiv($x, 0x10000);
  400.         ($y, $yr) = &bdiv($y, 0x10000);
  401.         $r = &badd(&bmul(int $xr | $yr, $m), $r);
  402.         $m = &bmul($m, 0x10000);
  403.     }
  404.     $r;
  405.     }
  406. }
  407.  
  408. # compute x ^ y
  409. sub bxor { #(num_str, num_str) return num_str
  410.     local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
  411.     if ($x eq 'NaN' || $y eq 'NaN') {
  412.     'NaN';
  413.     } else {
  414.     while ($x ne '+0' || $y ne '+0') {
  415.         ($x, $xr) = &bdiv($x, 0x10000);
  416.         ($y, $yr) = &bdiv($y, 0x10000);
  417.         $r = &badd(&bmul(int $xr ^ $yr, $m), $r);
  418.         $m = &bmul($m, 0x10000);
  419.     }
  420.     $r;
  421.     }
  422. }
  423.  
  424. # represent ~x as twos-complement number
  425. sub bnot { #(num_str) return num_str
  426.     &bsub(-1,$_[$[]);
  427. }
  428.  
  429. 1;
  430. __END__
  431.  
  432. =head1 NAME
  433.  
  434. Math::BigInt - Arbitrary size integer math package
  435.  
  436. =head1 SYNOPSIS
  437.  
  438.   use Math::BigInt;
  439.   $i = Math::BigInt->new($string);
  440.  
  441.   $i->bneg return BINT               negation
  442.   $i->babs return BINT               absolute value
  443.   $i->bcmp(BINT) return CODE         compare numbers (undef,<0,=0,>0)
  444.   $i->badd(BINT) return BINT         addition
  445.   $i->bsub(BINT) return BINT         subtraction
  446.   $i->bmul(BINT) return BINT         multiplication
  447.   $i->bdiv(BINT) return (BINT,BINT)  division (quo,rem) just quo if scalar
  448.   $i->bmod(BINT) return BINT         modulus
  449.   $i->bgcd(BINT) return BINT         greatest common divisor
  450.   $i->bnorm return BINT              normalization
  451.   $i->blsft(BINT) return BINT        left shift
  452.   $i->brsft(BINT) return (BINT,BINT) right shift (quo,rem) just quo if scalar
  453.   $i->band(BINT) return BINT         bit-wise and
  454.   $i->bior(BINT) return BINT         bit-wise inclusive or
  455.   $i->bxor(BINT) return BINT         bit-wise exclusive or
  456.   $i->bnot return BINT               bit-wise not
  457.  
  458. =head1 DESCRIPTION
  459.  
  460. All basic math operations are overloaded if you declare your big
  461. integers as
  462.  
  463.   $i = new Math::BigInt '123 456 789 123 456 789';
  464.  
  465.  
  466. =over 2
  467.  
  468. =item Canonical notation
  469.  
  470. Big integer value are strings of the form C</^[+-]\d+$/> with leading
  471. zeros suppressed.
  472.  
  473. =item Input
  474.  
  475. Input values to these routines may be strings of the form
  476. C</^\s*[+-]?[\d\s]+$/>.
  477.  
  478. =item Output
  479.  
  480. Output values always always in canonical form
  481.  
  482. =back
  483.  
  484. Actual math is done in an internal format consisting of an array
  485. whose first element is the sign (/^[+-]$/) and whose remaining 
  486. elements are base 100000 digits with the least significant digit first.
  487. The string 'NaN' is used to represent the result when input arguments 
  488. are not numbers, as well as the result of dividing by zero.
  489.  
  490. =head1 EXAMPLES
  491.  
  492.    '+0'                            canonical zero value
  493.    '   -123 123 123'               canonical value '-123123123'
  494.    '1 23 456 7890'                 canonical value '+1234567890'
  495.  
  496.  
  497. =head1 Autocreating constants
  498.  
  499. After C<use Math::BigInt ':constant'> all the integer decimal constants
  500. in the given scope are converted to C<Math::BigInt>.  This conversion
  501. happens at compile time.
  502.  
  503. In particular
  504.  
  505.   perl -MMath::BigInt=:constant -e 'print 2**100'
  506.  
  507. print the integer value of C<2**100>.  Note that without conversion of 
  508. constants the expression 2**100 will be calculated as floating point number.
  509.  
  510. =head1 BUGS
  511.  
  512. The current version of this module is a preliminary version of the
  513. real thing that is currently (as of perl5.002) under development.
  514.  
  515. =head1 AUTHOR
  516.  
  517. Mark Biggar, overloaded interface by Ilya Zakharevich.
  518.  
  519. =cut
  520.